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1.  BACKGROUND 


A.  HISTORICAL  DEVELOPMENT  OF  GRAY  EOS 

The  original  GRAY  three-phase  equation  of  state  (EOS)  for  metals  was  developed 
by  Royce1  at  the  Lawrence  Livermore  Laboratory(LLL).  It  included  a  Gruneisen 
description  of  the  solid  referenced  to  the  solid  Hugoniot2,  and  a  scaling  law  EOS  for 
liquid  developed  by  Grover2.  Material  in  the  liquid-vapor  region  was  described  by  a 
classical  hard-sphere  EOS  due  to  Young  and  Alder4  which  included  a  van  der  Waals 
attractive  term. 

The  solid-liquid  model  assumed:  the  entropy  of  the  melting  is  independent  of 
pressure;  the  temperature  dependence  of  the  specific  heat  in  the  liquid  is  a  universal 
curve  scaled  on  the  melting  temperature;  and  the  pressure  dependence  of  the  melting 
temperature  is  given  by  a  modified  Lindemann  law. 

The  Grover  scaling  law  and  the  Young-Alder  model  were  analytically  joined  at  a 
volume  in  the  range  1.3  to  1.5  times  the  normal  volume  by  adding  correction  terms  to 
the  Young-Alder  model  at  the  low  density  side  of  the  join  volume. 

Young5  subsequently  introduced  a  modification  of  the  original  GRAY  EOS  to 
permit  the  calculation  of  the  pressure  in  the  mixed-phase  liquid-vapor  region;  however 
the  basic  physical  models  and  join  procedure  were  not  altered. 

A  substantially  modified  GRAY  EOS,  BRLGRAY,  was  developed  in  1976  under 
Ballistic  Research  Laboratory  (BRL)  sponsorship  by  Wray6  at  Systems,  Science  and 
Software.  Modifications  included  a  continuous  model  for  the  liquid  phase,  a  soft- 
sphere  model  for  the  vapor  region,  and  an  improved  join  procedure  for  the  liquid  and 
vapor  regions. 

1.  E.B.  Royce.  "GRAY,  A  Three-Phase  EqUdtion  qf State  Jbr-MetalA  Lawrence  Livermore.  . 
Laboratory.  UCRL-51121,  Sep  1971. 

2.  R.  N.  Keeler  and  E.  B.  Royce,  " Shock  Waves  in  Condensed  Media ",  Lawrence  Liver¬ 
more  Laboratory,  UCRL-71846,  1969. 

3.  R.  Grover,  "  Liquid  Metal  Equation  of  State  Based  on  Scaling,  J.  Chem  Phys,  55_, 
3435.  1971. 

4.  D.  A.  Young  and  B.  J.  Alder,  " Critical  Point  qf  Metals  from  the  van  der  Waals  Moder , 

Phys  Rev  A3,  364,  1971.  '  ~ 

J.  D.  A.  Young.  "Modification  of  the  GRAY  Equation  of  State  in  the  Liquid-Vapor  Region" 
Lawrence  Livermore  Laboratory,  U  CRL-5 1575 Apr~T5 , 19T4. 

6.  William  O.  Wray  and  Robert  A.  Cecil.  "Modified  GRAY:  An  Improved  Three-Phase 
Equation  qf  State  for  Metal  f .  BRL  CR  No.  299,  April  1976.  (AD  #A025260) 
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The  fluid  (liquid  and  gas)  models  share  a  single  continuous  representation  of  the 
specific  heat.  This  eliminates  the  artificial  hot  liquid  region  of  the  original  model. 
The  soft  sphere  representation  in  the  vapor  region  yields  a  model  which  is  physically 
realistic  at  higher  temperatures  and  more  compatible  with  the  Grover  model. 

B.  APPLICATIONS 

Multi-phase  EOS’s  are  required  for  the  numerical  modelling  of  the  response  of 
materials  to  dynamic  loading  under  which  solid  phase  changes,  melt  transitions,  and 
vapor  states  occur.  However,  equations  of  state  currently  in  use  for  such  problems 
include  the  Tillotson  EOS7  which  is  used  in  the  hydrocodc  HELP758  ,  and  curve  fits 
such  as  that  used  in  the  HULL9  code.  These  EOS’s  do  not  take  into  account  phase 
changes  which  may  be  important  in  calculating  material  response  to  dynamic  loading. 

HELP75  is  the  most  recent  documented  version  of  the  HELP  code.  However,  it 
does  not  calculate  internal  energy  accurately  and  therefore  cannot  properly  predict 
thermodynamic  state.  HULL  is  an  Eulerian  code  which  uses  a  curve  fit  to  EOS  data, 
but  it  uses  the  same  type  of  energy  calculation  as  does  HELP75. 

BRLHELP10  is  a  new  improved  version  of  HELP.  It  has  a  modification  due  to 
Schmitt1 1  which  enables  it  to  calculate  internal  energies  very  accurately,  has  improved 
input  formatting12  and  in  its  latest  version  has  the  BRLGRAY  EOS  for  metals. 

Janet  Lacetera13  has  done  parameter  studies,  using  the  HELP  code,  for  a  variety 


7.  J.  H.  Tillotson,  " Metal  Equations  of  State  for  Hypervelocity  Impact ",  General  Atomic 
Report  GA-3216.  July  1962. 

8.  Laura  J.  Hageman,  et  al,  "HELP,  a  Multi-Material  Eulerian  Program  for  Compressible 
•  .  .  Fluid  and  Elastic-Plqstic  glows  in  Two  Space  Dimensions  and  Time”,  Systems,  Science,  and 

Software  Report  TR-76-45-BK2,  Apr  1976. 

9.  Richard  E.  Durrett  and  David  A.  Matuska,  "The  HULL  Code,  A  Finite  Difference  Solu¬ 
tion  to  the  Equations  of  Continuum  Mechanics?,  Air  Force  Armament  Laboratory,  AFATL- 
TR-78-125,  Nov  1978. 

10.  J.  Lacetera,  J.  E.  Lacetera,  and  J.  A.  Schmitt,  "The  BRL  7600  Version  of  the  HELP 
Code”,  BRL  REPORT  ARBRL-TR-02209.  Jan  80.  (AD  #A082804) 

11.  J.  A.  Schmitt,  "An  Improved  Internal  Energy  Calculation  for  the  HELP  Code  and  Its 
Implications  to  Conical  Shaped  Charge  Simulation? ,  BRL  Report  ARBRL-TR-02168,  Jun 
79.  (AD  #A072785) 

12.  J.  Lacetera,  ” BRLHELP80 :  NAMELIST  INPUT",  BRL  IMR  No  677,  Mar  80. 

13.  Janet  Lacetera,  ” Study  of  Liner  Collapse  and  Jet  Formation  for  Various  Hemispherical 
Shaped-Charge  Designs”,  BRL  Draft  Report. 
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of  problems.  Agreement  with  experiment  for  these  calculations  has  been  quite  good 
and  indicates  that  the  HELP  code  provides  sufficiently  accurate  information  to  be  used 
as  a  design  optimization  tool  in  hydrocode  applications.  However,  the  EOS  imbedded 
in  the  code  does  not  adequately  represent  the  effects  of  phase  changes.  For 
applications  where  such  phenomena  are  important  a  more  sophisticated  EOS  is 
appropriate. 

We  have  made  comparisons  between  BRLGRAY,  Tillotson,  and  the  HULL 
EOS’s.  In  the  solid  state  the  agreement  between  these  EOS’s  and  with  experiment  is 
excellent.  Significantly,  BRLGRAY  uses  no  more  computer  time  than  do  the  others 
for  solid  state  calculations  and  outperforms  the  others  in  the  melt  transition  and  at 
higher  temperatures.  A  comparison  of  the  BRLGRAY  EOS  and  that  used  in  the 
HULL  code,  for  temperature  as  a  function  of  energy  for  copper,  is  shown  in  Figure  1 . 
The  phase  transition,  in  BRLGRAY,  begins  at  an  energy  of  4.62  ergs/g  and  is 
complete  at  an  energy  of  6.67  ergs/g.  This  heat  of  transition  of  2.05  ergs/g  is  very 
close  to  the  values  given  in  the  various  handbooks.  Figure  2  is  another  illustration  of 
the  BRLGRAY  melt-phase  transition:  ambient-pressure  temperature  profiles  are 
shown  for  the  copper  and  aluminum.  The  melt  temperatures  are  1356  K  for  copper 
and  932  K  for  aluminum.  The  melt  transition  is  simply  not  calculated  by  the  other 
EOS’s. 

2.  The  BRLGRAY  Equation  of  State* 

A.  Solid  Region 

This  portion  of  the  EOS  remains  unmodified.  Grover’s  scaling  law  Includes  a 
Gruneisen  description  of  the  solid  referenced  to  the  experimental  Hugoniot.  The 
volume  dependence  of  gamma  in  the  solid  is  linear.  Energy  and  pressure  in  this 
region  are  given  by 

P't2 

ES(T,V )  -  EC(V)  +  3R'  +  (1) 

and 

PS(T,V)  -  P0(V)  +  3ys(V)-^  4-  yye-~ -,  (2) 

where  R’  and  G’  are  the  atomic-weight  scaled  gas  constant  and  electronic  energy 
coefficient  respectively,  Po(V)  and  Eo(V)  are  reference  values  at  OP  K,  and  ys(V)  is 
the  volume-dependent  Gruneisen  coefficient  for  the  solid. 


*  A  list  of  symbols  for  the  EOS  is  given  in  Appendix  A, 
and  a  complete  list  of  the  EOS  equations  is  given  in  Appendix  B. 


TEMPERATURE  CK> 


The  atomic-weight  scaled  electronic  energy  coefficient  is  given  by 

rr2NekR'  ’ 

I  O 


where  the  Fermi  energy,  Ef  is  given  by 

Ef-26| 


2Ef  ’ 


N.  ]2/3 


where  W  is  atomic  weight,  and  Ne  is  number  of  free  electrons. 

The  complete  solid  EOS  is  described  by  Royce1.  The  experimental  Hugonlot  is 
assumed  to  have  the  standard  form: 

Us  -  C  +  SUp,  (5) 


PH(V) 


P0C2x 

(1-Sx)2’ 


Eh(V)  -  }pH(V)[V0-V]  +  EpH,  (7) 

where  x  is  a  compression  variable,  and  C  and  S  are  Hugonlot  parameters.  To 
reference  the  EOS  to  the  Hugoniot,  r(T,V)  is  defined  by  the  the  Gruneisen-  like 
equation 

P(T,V)  -  PC(V)  +  [e(T,V)  -  E0(V)],  (8) 

and  the  HuRoniot  in  tTie  "solid  is  written  *  •  •  ....  .  .. .  ....... 

PH-P0(V)  +  -^[EH(V>-Et,(V)],  (9) 

where  ys  is  the  volume-dependent  Gruneisen  parameter  in  the  solid. 


Equations  7-9  are  combined,  using  r(T=0,V)— ys(V)  to  give 

P(T,V)  -  P„(V)U  -  ^(V0-V)l  +  ^ - [ht.V)  - 

+  [r(T.v)-n(v»|(E(T,v)_Eo(v)] 

where  EoH  is  the  energy  at  the  foot  of  the  Hugoniot. 

P(T,V)  is  then  written  as  the  sum  of  two  parts 
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P(T,V)  -  Pg(V,E)  +  PC(T,V), 


where 


pg(v,e)-ph(v) 


1  - 


y,(v) 

2V 


(v„-v)  + 


y.(v) 


[e(t,v)  —  Ejjh] 


is  a  Gruneisen-like  equation  of  state  for  the  solid,  and 

PC(T,V)-  >s(V)  |  [e(T,V)  -  Et)(V )] 


(11) 

(12) 

(13) 


is  a  correction  term. 

Since  the  EOS  is  referenced  to  the  experimental  Hugoniot,  an  explicit 
expression  for  Po(V)  is  not  needed.  However,  the  calculation  of  the  temperature,  as 
a  function  of  energy,  which  is  needed  for  the  calculation  of  Pc  requires  a  knowledge 
of  E0(V).  This  is  given  by  Royce1  to  be 


C2  —  x2 
2(1— Sx) 


c  .  s  ,  s2  . 
s  +  T  +  T 


y0 


)x2 


(14) 


which  he  estimates  to  be  correct  to  +10  percent  for  compressions  up  to  2. 


B.  Melt  and  Liquid  Regions 


The  Grover  liquid-metal  EOS  is  based  on  the  observation  that  the  entropy  of 
melting  can  be  scaled  on  the  atomic  weight  using  a  single  constant  (2.32  cal/g-atom 
deg)  for  all  metals  and  that  the  liquid  specific  heat  has  a  universal  dependence  on  the 
scaled  .terqperatyreCT/Tm)-  These  observations  are  supported  by  experimental  data. 
Finally,  the  Lindemann  law  is  used  to  establish  the  volume  dependence  of  the  melting 
temperature. 


In  BRLGRAY  the  fluid  specific  heat  has  a  universal  hyperbolic  dependence  on 
the  scaled  temperature: 


Cv 


3R 

2 


T 

-r 

1  + 

a— — H  1 

T 

a*=0.1,  T>Tm 


(15) 


This  is  illustrated  in  Figure  3,  along  with  the  Grover  original  linear  specific  heat  for 
the  liquid.  The  Wray  form  limits  asymptotically  to  the  ideal  gas  value. 


The  entropy  in  the  liquid  region  is  obtained  from  the  thermodynamic  relation 

as  cv 

3T  |y  “  T  ’ 

Using  Equation  5  and  integrating  Equation  6  with  respect  to  temperature  gives, 
following  Wray’s  notation. 
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Sj  —  3RlnT  —  -^-ln  a~~  +  1  +fXV) 


Similarly  Tor  the  solid, 


Sg  -  3RlnT  +  h(V), 


where  f  and  h  are  undetermined  functions  of  volume.  The  entropy  of  melting  is  then 
defined  as 

AS  -  fTV)  -  h(V), 

where  AS  is  determined  from  experimental  data2. 

The  entropy  in  the  liquid  is  then 

Sj  ==  Ss  ~~ln  a~~  -1"  1  +  AS. 


In  the  two-phase  melt  region. 


Sm  ■=*  vS!  +  (1— v)Ss 


where 


T  -  TS(V) 


,  and  T„  — 


Ts  +  T, 


"  T|(V)— TS(V)  m  2 

Since,  in  the  melt  region,  T  ~  Tm,  the  meit  entropy  is  approximated  by 


Sm  “  Ss  +  vlAS  —  .1 43R|. 


The  Helmholtz  free  energy  A  is  obtained  from  the  thermodynamic  relation 


S«  — 


flA 

ST  L’ 


which  may  be  integrated  to  give 

A-E^V)-  /SdT, 

where  Eo(V)  is  the  cold  compression  energy.  This  yields,  after  substitution  and 
integration, 

ATO-As=*»'2(AT/2)(AS-.143R),  07) 

where  AT  -  T,(V)  -  TS(V)  and  from  A  -  E  -  TS  it  follows  that 

Em^Es  +  v[T-  »(AT/2)](AS  -  .143R),  08) 

where  Es'=As+TSs  is  the  extrapolated  energy  in  the  solid. 

The  pressure  in  the  melt  region  is  obtained  by  applying  the  relation 


p  -  - 


8A 


av 


to  Equation  17,  and  using  the  approximations 

dUn.Q  dT»^dTm 
dV  ’  dV  dV 

along  with  the  volume  dependence  of  the  melt  temperature,  given  by  the  Lindemann 
law, 


X  - 


dlnV 


to  get 


Pm  -  Ps  +  -  .143R). 


(19) 


The  same  procedure,  starting  with  the  entropy  in  the  liquid  region  yields  similar 
equations  for  the  liquid  region.  These  are: 


c  c  3R, 
si  -  Ss  -  —In 


«T 


*  m 


+  1 


+  AS; 


(20) 


A,  -  As  -  -^p(AS  -  .143R)  -  (T  -  Tm) 
3RT, 


AS  4- 


3R 


1  m 


2a 


[•£■+1 

In 

•£+i 

—  (a  +  l)ln(a  +  1)  ; 

lT» 

T 

(21) 


3RT. 


2a 


T  .  . 

In 

a  — - h  1 

T 

E,-Es  +  TniAS  +  .143R^ 
—(a  +  l)ln(a  +  1)  — 


T  i 

T"1 


;and 


(22) 


Pi- 


3RX  Tm 
2a  V 


In 

Lx+il 

T 

1  m 

—  (a+l)ln(a+l)  — 


All  of  these  equations  are  listed  for  reference  in  Appendix  B. 


(23) 
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C.  The  Vapor  Region 


The  classical  van  der  Waals  equation  of  state  has  the  form 


P  -  RT 

v  V-Vh 


where  a  is  the  coefficient  of  attractive  potential,  and  Vb  is  the  hard  sphere  excluded 
volume  which  represents  the  repulsive  potential  by  increasing  the  pressure  at  small 
volumes.  The  Young-Alder  modified  form  of  the  van  der  Waals  EOS  has  the  form 


(l-r,)3 


where  tj  *=  — . 

The  first  term  in  Equation  24  is  an  algebraic  approximation14  to  the  hard  sphere 

pressure,  which  was  developed  from  an  analysis15  of  the  reduced  virial  series,  and  has 

been  shown  to  describe  rigid-sphere  behaviour  better  than  any  other  existing  form. 

An  ideal  gas  caloric  equation  of  state: 

„  3RT  _a_ 

Ev-—  v 

completed  the  model. 

Wray6  noted  that  the  Gruneisen  coefficient. 


r  =  2 
G  3 


(l— n)3 


in  the  Young-Alder  model  yielded  an  unreasonably  high  value  at  the  join  volume 
between  the  Grover  and  Young-Alder  models,  and  being  independent  of  temperature 
did  not  limit  properly  as  temperature  approached  infinity.  His  soft-sphere  model 
allows  Vb  to  vary  as  a  function  of  temperature,  such  that 


'  Vi,(T), 


lim  Vb 

T— >°o 


T— oo  dT 


thus  yielding  the  correct  ideal  gas,  high  temperature  limit  of  G-2/3. 


14.  N.  F.  Carnahan  and  K.  E.  Starling,  Journal  of  Chemical  Physics,  51,  635,  1969. 

15.  B.  J.  Aider  and  T.  E.  Wainwright,  Journal  qf  Chemical  Physics.  33,  1439,  1960. 
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The  generalized  Rowlinson  method16  was  used  to  relate  the  temperature 
dependence  of  the  sphere  diameter  to  realistic  fluid  behaviour,  so  the  soft  sphere 
diameter  is  calculated  according  to 


d  - 


where  u(r)  is  the  repulsive  potential.  Since  Grover  showed  his  liquid  scaling  law  to  be 
in  close  agreement  with  computer  experiments  in  the  liquid  phase  for  particles  having 
an  inverse-twelve  potential,  the  following  potential  is  used, 


u(r)  - 


12 


where  e  is  the  Fermi  energy,  and  o-  is  the  hard-sphere  diameter, 
eventually, 


d 


1.0556 


ter 


12 


kT 


1/12 


This  yields. 


for  the  soft-sphere  diameter. 

Using  the  appropriate  thermodynamic  relations,  the  thermal  energy  E  and  the 
Helmholtz  free  energy  A  can  be  derived  as 


3RT  2RT2  (2 — 17)  dVb  a 

^“2  V  (1 — »j)3  dT  V  ^ 

and 


(25) 


A. 


-  lnV 

2(1— rj)3 


lnT-^  +  E, 


The  complete  vapor  EOS  for  BRLGRAY  is  listed  in  Appendix  B. 


(26) 


D.  Join  Region 

In  Wray’s  join  region,  between  the  liquid  and  vapor  states,  a  trigonometric  mixing 
function  was  applied  to  the  Helmholtz  free  energy  calculated  in  each  of  the  two  states. 
All  other  functions  of  interest  were  then  obtained  through  the  application  of  the 


16.  H.  L.  Anderson,  J.  D.  Weeks  and  D.  Chandler,  " Relationship  between  the  Hard- 
Sphere  Fluid  and  Fluids  with  Realistic  Repulsive  Forces”,  Physical  Review,  Vol.  4,  No.  4, 
October  1971. 


> 
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appropriate  thermodynamic  relations.  However,  this  method  is  replaced  in 
BRLGRAY  with  a  simple  volume  weighted  interpolation  scheme  which  gives  better 
join  region  pressures  by  avoiding  less  accurately  calculated  functions  in  the  join  region. 

3.  EQUATION  OF  STATE  INPUT  PARAMETERS 

The  EOS  input  parameters  are  stored  in  a  two-dimensional  array, 
ESTCON(I,N), where  1  is  the  parameter  index  and  N  is  the  material  index.  This  array 
is  defined  in  Table  1.  Table  2  lists  the  ESTCON  data  for  the  metals  Al,  Cu,  and  Fe. 
It  should  be  noted  that  V,  and  V2  in  the  ESTCON  array  (1—14  and  1—15)  are 
defined  in  the  BRLGRAY  EOS  and  represent  the  volume  limits  of  the  new  join 
region.  They  replace  the  single  join  volume,  Vj(I— 25), of  the  original  GRAY  EOS. 

The  ESTCON  array  is  maintained  in  data  statements  in  a  subroutine  named 
ESTKON.  Note  that  Ae  has  the  constant  value  .6667  for  all  metals  and  AS  has  the 
constant  value  9.637E-05  for  all  metals. 


TABLE  I 

ESTCON  PARAMETERS 

Index 

Symbol 

Definition  "1 

1 

Po 

Normal  density  Kg/m3 

2-5 

Unused 

6 

cQ 

Bulk  Sound  Speed  (cm/s) 

7 

AMU 

Shear  Modulus  (d/cm2) 

8 

Y0 

Yield  Strength  (d/cm2) 

9 

<T 

Spall  Strength  (d/cm2) 

10 

E* 

Sublimation  Energy  (erg/g) 

11 

ELv 

Incipient  Vaporization  Energy  (erg/g) 

12 

Elm 

Complete  Melt  Energy  (erg/g) 

13 

Esm 

Incipient  Melt  Energy  (erg/g) 

14 

V, 

Liquid-Join  Boundary  (m3/Kg) 

15 

V2 

Join-Vapor  Boundary  (m3/Kg) 

16 

s 

Hugoniot  Parameter:  (cm/^sec) 

17 

*0 

Lattice  gamma:  A(V)— kQ— ax 

18 

a 

Defined  in  ESTCON(17,N) 

19 

*e 

Electronic  gamma 

20 

Be 

Electronic  energy  coefficient 
(mbar  cm3/mole  deg2) 

21 

tmo 

Melt  temperature  (°K)  at  V— VQ 

22 

EoH 

Energy  at  V-Vo>T-300PK, 

P—0, (mbar  cm2/g) 

23 

Eoo 

Energy  difference  between  E^  and 
energy  at  V«*Vo,T-0(mbar  cm3/g) 
Default  value:  E^— 

-3000*8.1 34*10~s+l  50g,,)/W 

24 

w 

Atomic  weight  (g/mole) 

25 

Unused 

26 

Vb 

Vapor  exclusion  volume  (m3/Kg) 

27 

aY 

Coefficient  of  attractive  potential 
mbar(cm3/mole)2 

28 

Tc 

Critical  temperature  (°K) 

29 

AS 

Entropy  of  melting  mbarcm3/mole°) 

30 

EOSFLG 

Equation  of  state  flag 

21 


ES’fe*  ■ 


TABLE  2 


ESTCON  DATA 


Symbol 

A1 

Cu 

Fe 

Po 

2.7 

8.94 

7.87 

C0 

5.44E+05 

3.94E+05 

5.13E-M55 

E* 

1.18E+11 

S.37E+10 

6.80E+10 

ELv 

2.6SE+10 

4.79E+10 

Elm 

1.06E+I0 

6.67E+09 

1.49E+10 

Esm 

6.62E+09 

4.62E+09 

1.16E+10 

v, 

.40 

.12 

.14 

V2 

.58 

.16 

.18 

s 

1.34 

1.49 

1.40 

2.18 

1.97 

1.80 

a 

1.70 

1.50 

1.50 

*e 

.6667 

.6667 

.6667 

8e 

8.70E-09 

4.90E-09 

9.60E-09 

Tmo 

1220 

1790 

2385 

w 

26.98 

63.54 

55.85 

Vb 

.190 

.056 

.066 

aY 

47.0 

35.0 

31.0 

T 

*c 

10633 

11407 

9753 

4.  THE  DATGEN  PROGRAM 

In  order  to  use  the  BRLGRAY  EOS  in  the  liquid-vapor,  mixed-phase  region  it  is 
necessary  to  calculate  the  four  critical  constants: 

Vc,  Tc,  and  PC(VC,TC)  and  Et(Vc,Tc), 

and  the  mixed-phase  boundary.  We  calculate  these  using  a  program  which  is  based  on 
the  main  subroutine  given  by  Young5 . 

Initially  the  critical  volume  is  estimated  to  be  4.3S  times  the  boiling  point  volume. 
An  iteration  process,  using  pressures  P(V,T),  generated  by  the  EOS,  then  locates  the 
critical  isotherm  and  the  volume  for  which  AP/AV  is  zero  on  the  critical  isotherm,  the 
critical  volume.  At  this  point  P— Pc  and  E^.  is  calculated  by  the  EOS.  The  code  then 
chooses  a  set  of  equally  spaced  temperatures  between  Tc  and  Tmo,  the  melting 


temperature  parameter,  to  draw  tie  lines  connecting  the  high  density  liquid  with  the 
low  density  vapor.  The  corresponding  pressures  and  equilibrium  liquid  and  vapor 
densities  are  computed  by  Maxwell  construction.  The  locus  of  these  tie-line  end 
points  then  forms  the  dome-shaped,  mixed-phase  boundary  in  PV  space.  The 
corresponding  energies  E(Vmin,T)  and  E(Vmax,T)  are  also  computed  for  each 
temperature  and  the  data  is  stored  in  a  permanent  file  for  access  by  the  EOS. 

The  associated  variable  names  are: 


TMP:temperature  on  a  tie  line  (°K) 

VMN:volume  on  liquid  end  of  tie  line  (  m3/Mg) 

VMX:volume  on  vapor  end  of  tie  line  (  m3/Mg) 

EVMN:energy  at  VMN  and  TMP  (  erg/g  ) 

EVMX:energy  at  VMX  and  TMP  (  erg/g  ) 

PRES:pressure  on  a  tie  line  (  d/cm2  ) 


where  d/cm2  are  working  units  in  the  code,  and  are  equivalent  to  10  1  Pascals. 

The  DATGEN  program  also  has  an  isotherm  driver  subroutine  which  generates 
isotherms  in  P-V  space  using  the  EOS  along  with  subroutine  ET5.  Subroutine  ET 
calculates  the  specific  internal  energy  consistent  with  the  driver  values  of  specific 
volume  and  temperature.  The  BRLGRAY  EOS  is  then  called  to  calculate  the 
pressure.  This  program  is  used  to  determine  the  volume  limits  of  the  new  join 
region  in  BRLGRAY.  Wray6  performed  such  calculations  for  three  metals:  Al,  Fe, 
and  U.  We  have  made  calculati&ns' for'  seventeen  materials  using  slightly  different, 
procedures,  and  in  this  report  have  presented  data  for  Al,  Fe,  and  Cu.  These 
calculations  are  described  in  the  next  section  of  this  report. 


5.  MIXED-PHASE  DATA  GENERATION 


Figure  4  shows  a  P-V  plot  of  the  regions  defined  in  BRLGRAY.  Pressure  and 
temperature  from  the  liquid  and  vapor  phases  are  volume  weighted  in  the  join  region. 

Figure  5  shows  a  schematic  diagram  of  representative  isotherms  in  the  join 
region.  The  values  of  V,  and  V2  are  chosen  such  that  pressure  maxima  and 
minima  do  not  occur  along  an  isotherm  in  the  join  region.  Such  constructs  would 
make  it  impossible  for  the  DATGEN  program  to  find  the  critical  isotherm  in 
subsequent  mixed-phase  boundary  calculations. 

Wray  used  a  trial  and  error  procedure  in  which  isotherms  were  generated  for 
various  values  of  V,  and  V2  .  We  used  a  similar  isotherm  driver  but  ran  separate 
cases  in  which  the  EOS  was  constrained  to  either  the  solid-liquid  or  the  vapor  model. 
Analysis  of  these  pressure  isotherms  then  allowed  a  good  initial  estimate  of  Vj  and 

v2- 


The  isotherms  generated  by  the  solid-liquid  model  were  well  behaved  up  to  a 
volume  at  which  the  model  failed  catastrophically.  These  isotherms  gave  a  good 
indication  of  the  maximum  value  which  V,  could  take.  Comparison  of  these 
isotherms  with  isotherms  generated  by  the  vapor  model  then  allowed  an  estimate  of 
the  best  value  of  V2  . 

Isotherms  generated  using  these  values  of  Vt  and  V2  generally  gave  positive 
pressures  in  the  join  region,  without  maxima  or  minima,  for  the  temperatures  of 
interest. 

Alternatively,  once  these  parameters  have  been  generated  for  a  number  of  metals 
one  may  estimate  Vt  and  V2  for  a  new  metal  using  the  following  scaling  rule: 

PoV,-l.l«Cj 

PDV2  —  1.4  —  C2 

where  the  Cj  are  average  values  of  p0Vj.  The  V(  estimated  in  this  way  will  need  very 
little  adjustment  after  observation  of  the  isotherms  in  the  join  region  using  the 
complete  EOS. 

The  driver  routine  was  programmed  to  generate  isotherms  which  ranged  from 
20,000°K  down  to  2,000°K  in  decrements  of  2,000°K,  and  to  lower  temperatures  in 
smaller  steps.  Volumes  ranged  from  .01  m3/Mg  to  beyond  the  mixed-phase  region  in 
equal  increments  of  In  V. 

Aluminum  isotherms  were  generated  starting  with  a  specific  volume  of  .01  and 
incremented  until  the  calculation  failed.  Inspection  of  the  code  output  then 
determined  maximum  values  of  Vj  for  each  isotherm.  These  values,  shown  in 
Table  3,  ranged  from  .67  m3/Mg  at  20,000°K  to  .37  m3/Mg  at  300°K.  In 
anticipation  of  being  able  to  extend  the  validity  of  the  model  to  below  1,000°K,  at 
these  volumes,  we  chose  Vj  to  be  40  m3/Mg. 

Vapor-model  aluminum  isotherms  were  then  generated  starting  with  specific 
volumes  equal  to  .355  m3/Mg  and  extending  beyond  the  expected  value  of  V2  . 
Comparison  of  these  isotherms  with  those  generated  with  the  solid-liquid  model 
showed  a  crossover  point  for  each  isotherm  which  was  taken  as  the  minimum  value  of 
V2  in  each  case.  These  values,  shown  in  Table  4,  ranged  from  .58  m3/Mg  at 
20,000°K  to  .54  m3/Mg  at  6,000°K.  Thus  .58  m3/Mg  was  interpreted  as  an  upper 
bound  on  V2  ,  ensuring  that  crossover  would  occur,  within  the  join  region,  for  all 
temperatures  below  20,000°K. 

Finally,  a  set  of  isotherms  were  generated  in  the  code’s  normal  mode  (in  which  it 
determines  the  thermodynamic  phase  and  chooses  the  corresponding  model  of  the 
EOS)  using  the  values  V,-.40  m3/Mg  and  V2-.58  m3/Mg.  These  isotherms  gave 
positive  pressures  in  the  join  region  without  maxima  or  minima  for  temperatures 
down  to  1,000°K.  A  set  of  isotherms  using  V,  -  .41  m3/Mg  showed  a  minimum  in 
the  join  region  at  4,000°K,  and  a  set  using  .37  m3/Mg  showing  minima  in  the 
I,000°K  isotherm. 
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TABLE  3 

_ ^lmax 

for  Aluminum  Isotherms  Using  Solid-Liquid  Model 

Temperature 

^lmax 

(  nr/Mg  ) 

Phase 

(K) 

20,000 

.67 

15,000 

.66 

Liquid 

10,000 

.62 

8,000 

.61 

Melt 

6,000 

.58 

4,000 

.54 

2,000 

.45 

Solid/Melt/Liquid 

1,000 

.41 

500 

.37 

300 

.37 

v, 

TABLE  4 

min  For  Aluminum  Isotherms  Using  Vapor  Model 

Temperature 

V2min  Pha* 

(K) 

(  m3/Mg  ) 

I 

.58 

.56  Vapor 

.55 

.54 

Vj  and  V2  were  determined,  using  the  same  approach,  for  the  other  metals. 
These  are  shown  in  Table  5  along  with  the  values  for  Al. 

Once  the  volume  limits  of  the  join  region  were  determined  for  each  of  the  metals, 
the  DATGEN  program  was  used  to  calculate  the  critical  point,  and  subsequently  the 
mixed-phase  boundary.  Critical  point  data  for  the  five  metals  is  shown  in  Table  6. 
Finally  the  DATGEN  program  is  used  to  generate  the  mixed-phase  boundaries. 
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TABLE  5 

JOIN  VOLUMES 


Metal 

Vi 

v2 

A! 

.40 

.58 

Cu 

.12 

.16 

Fe 

.14 

.18 

TABLE  6 

CRITICAL  POINTS 


Metal 

Vc 

(  m3/Mg  ) 

Tc 

(°K) 

Pc 

(  lO8*?  ) 

Ec 

(.lmJ/g) 

A1 

1.460 

7180 

5.34E-03 

1.25E+09 

Cu 

.413 

8880 

9.88E-03 

6.34E+10 

Fe 

.498 

9910 

1.04E-02 

7.24E+10 

TABLE  7 

ALUMINUM  MIXED-PHASE  BOUNDARY 

T 

*c 

7178 

vc 

1.46 

Pc 

5.34E+09 

Ec 

1.25E+1 1 

TMP 

VMN 

VMX 

EVMN 

EVMX 

PRES 

6.64E+03 

8.93E-01 

3.03E+00 

1.07E+11 

1.35E+11 

3.83E+11 

6.09E+03 

7.63E-01 

4.76E+00 

9.87E+10 

1.37E+11 

3.83E+09 

5.55E+03 

6.85E-01 

7.69E+00 

9.16E+10 

1.38E+11 

1.68E+09 

5.01  E+03 

6.33E-0I 

1.25F.+01 

8.55E+10 

1.37E+11 

9.98E+08 

4.47E+03 

5.95E-01 

2.00E+0I 

8.01E+10 

1.36E+1 1 

5.31E+08 

3.93E+03 

5.68E01 

3.33E+01 

7.50E+10 

1.35E+11 

2.40E+08 

3.39E+03 

5.56E-01 

1.11E+02 

7.01E+10 

1.33E+I1 

8.86E+07 

2.35E+03 

5.46E01 

2.50E+02 

6.52E+10 

1.31E+11 

3.21E+07 

2.30E+03 

5.38E01 

1.11  E+03 

5.8 1  E-h  10 

1.28E+11 

5.70E+06 

1.76E+03 

5.10E-01 

1.43E+04 

4.14E+10 

1.26E+11 

3.64E+05 

1.22E+03 

4.22E-01 

3.33E+06 

1.40E+10 

1.24E+11 

9.88E+02 

TABLE  8 

COPPER  MIXED-PHASE  BOUNDARY 

84^2 

vc 

.435 

Pc 

8.91  E+09 

Ec 

6.33E+11 

TMP 

VMN 

VMX 

EVMN 

EVMX 

PRES 

7.83E+03 

2.69E-01 

8.93E-01 

5.49E+10 

6.85E+10 

6.50E+09 

7.22F.+03 

2.31  E-01 

1.39E+00 

5.08E+10 

6.96E+10 

4.S5E+09 

6.62E+03 

2.09E-01 

2.13E+00 

4.74E+10 

6.99E+10 

3.02E+09 

6.01  E+03 

1.93E-01 

3.33E+00 

4.45E+10 

6.98E+10 

1.87E+09 

5.41  E+03 

1.81  E-01 

5.55E+00 

4.18E+10 

6.93E+10 

1.05E+09 

4.81  E+03 

1.72E-01 

1.00E+01 

3.94E+10 

6.86E+10 

5.23E+08 

4.20E+03 

1.65E-01 

2.00E+01 

3.73E+10 

6.78E+I0 

2.I7E+08 

3.60E+03 

1.59E-01 

5.00E+01 

3.52E+10 

6.68E+10 

6.93E+07 

2.98E+03 

1  55E-01 

2.00E+02 

2.94F.+  10 

6.57E+10 

1.74E+07 

2.39E+03 

1.49E-01 

1.67E+03 

2.28E+10 

6.46E+10 

1.83E+06 

1.79E+03 

1.38E-01 

5.00E+04 

1.33E+10 

6.34E+10 

3.50E+04 
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TABLE  9 


IRON  MIXED-PHASE  BOUNDARY 


T 

*c 

9913 

vc 

.498 

Pc 

1.04E+10 

Ec 

7.23E+10 

TMP 

VMN 

VMX 

EVMN 

EVMX 

PRES 

9.23E+03 

3.12E-01 

1.01E+00 

6.15E+10 

7.94E+10 

7.69E+09 

8.54E+03 

2.67E-01 

1.54E+00 

5.60E+10 

8.08E+10 

5.46E+09 

7.86E+03 

2.41E-01 

2.33E+00 

5.16E+10 

8.13E+10 

3.70E+09 

7.18E+03 

2.23E-01 

3.57E+00 

4.78E+10 

8.11E+10 

2.35E+09 

6.49E+03 

2.09E-01 

5.88E+00 

4.43E-M0 

8.06E+10 

1.38E+09 

5.81E+03 

1.99E-01 

1.00E+01 

4.11E+10 

7.98E+10 

7.20E+08 

5.12E+03 

1.90E-01 

2.00E+01 

3.82E+10 

7.87E+10 

3.24E+08 

4.44E+03 

1.83E-01 

5.00E+01 

3.36E+10 

7.75E+10 

1.14E+08 

3.75E+03 

1.73E-01 

1.67E+02 

2.63E+10 

7.61E+10 

3.14E+07 

3.07E+03 

1 .49E-01 

1.11E+03 

1.76E+10 

7.47E+10 

3.83E+06 

2.39E+03 

1.43E-01 

2.00E+04 

1.36E+10 

7.32E+10 

1.43E+05 

6.  Conclusion  j; 

I  ’■ 

Equatlon-of-state  input  parameters  have  been  generated  for  seventeen  different  | 

materials.  Data  has  been  shown  here  for  Al,  Cu,  and  Fe.  The  EOS  has  been 
installed  in  the  BRLHELP  code  and  comparisons  have  been  made  with  the  j 

Tillotson  EOS  and  the  temperature  fit  in  the  HULL  code.  In  the  solid  state  the 
agreement  between  BRLGRAY,  Tillotson,  the  HULL  EOS’s,  and  with  j< 

experimental  data  is  excellent.  Significantly,  BRLGRAY  uses  no  more  computer 
time  for  the  solid  state  than  do  the  other  EOS’s,  while  outperforming  them  in 
calculating  phase  change  phenomena. 

The  BRLGRAY  EOS  is  being  made  available  tor  use  by  hydrocode  users  at  the 
BRL  and  at  the  ARC  computers  at  Huntsville,  as  part  of  the  BRLEOS  code 

system17.  j 
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APPENDIX  A 


All  symbols  used  in  this  report  are  defined  where  they  first  appear  in  the  text 
However,  selected  symbols  are  listed  here  lor  convenience. 


(A-l) 

R’  -  R/W 

atomic  weight  scaled  gas  constant 

(A-2) 

G’  -  ge  /W 

atomic  weight  scaled  electronic  energy  coefficient 

(A-3) 

ys(v) 

volume  dependent  Gruneisen  coefficient  for  solid 

(A-4) 

ye 

electronic  gamma  (normally  2/3) 

(A-5) 

w 

atomic  weight 

(A-6) 

Ne 

number  of  free  electrons 

(A-7) 

us 

shock  velocity 

(A-8) 

Ur 

particle  velocity 

(A-9) 

C 

Us  intercept  at  zero  pressure 

(A-10) 

s 

slope  of  shock/particle  velocity  data 

(A-U) 

V-V 

x“  V 

compression  variable 

(A-12) 

Eoh 

Energy  at  V-Vn  ,T-300°K,P-0 

(A-13) 

Eoo 

Energy  difference  between  ILh  and  energy  at  V— VD  ,T—0 
Default-  -300(3*8.134*10-H-150gc)/W 

(A-14) 

a 

constant  in  hyperbolic  form  of  specific  heat 

(A-15) 

AS 

constant  entropy  of  melting 
constant— 9.62*10  5  Mbar  cnr/mole0 

(A-16) 

T 

1  m 

melting  temperature  (°K) 
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